Localization of shocks in driven diffusive systems without particle number 

conservation 



V. Popkov \ A. Rakos 1 , R. D. Willmann 1 , A. B. Kolomeisky 2 and G. M. Schiitz 1 
Institut fur Festkorperforschung, Forschungszentrum Julich - 52425 Jiilich, Germany and 
Department of Chemistry, Rice University, Houston, Texas 77005-1892 
(Dated: 6th February 2008) 

£SJ . We study the formation of localized shocks in one-dimensional driven diffusive systems with 

spacially homogeneous creation and annihilation of particles (Langmuir kinetics) .We show how to 
obtain hydrodynamic equations which describe the density profile in systems with uncorrelated 
steady state as well as in those exhibiting correlations. As a special example of the latter case 
the Katz-Lebowitz-Spohn model is considered. The existence of a localized double density shock 
is demonstrated for the first time in one-dimensional driven diffusive systems. This corresponds to 
phase separation into regimes of three distinct densities, separated by localized domain walls. Our 
analytical approach is supported by Monte-Carlo simulations. 
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I. INTRODUCTION 



One-dimensional driven diffusive systems proved to be a rewarding research topic in thepast years They were 
shown to exhibit boundary induced phase transitions |2| , spontaneous symmetry breaking |3|, |j] and phase separation 
0,0. Recently, attention turned to the case of systems without particle conservation in the bulk. In Ref.0the effect of 
a single detachment site in the bulk of an asymmetric simple exclusion process (ASEP) was studied. In Refs.0,0] the 
T \ interplay of the simplest one-dimensional driven model, the totally asymmetric exclusion process (TASEP) with local 
' O • absorption/desorption kinetics of single particles acting at all sites, termed 'Langmuir kinetics' (LK) was considered. 
These models were inspired by the dynamics of motor proteins |22|| which move along cytoskeletal filaments in a certain 
prefered direction while detachment and attachment can also occur between the cytoplasm and the filament, and, in 
a very different setting, by dynamics of limit orders in a stock exchange market. Being an equilibrium process, LK is 
well understood, while the combined process of TASEP and LK showed the new feature of a localized shock in the 
density profile of the stationary state |8| . 

The TASEP is defined on a one-dimensional lattice of size L. Each site can either be empty or occupied by one 
particle. In the bulk particles can hop from site i to site i + 1 with unit rate, provided the target site is empty. At 
| site 1 particles can enter the lattice from a reservoir with density p_ provided the site is empty. They can leave the 
system from site L into a reservoir of density p + with rate 1 — p+. Thus in the interior of the lattice the particle 
number is a conserved quanti ty. The phase diagram and steady states of the TASEP as a function of the boundary 
rates are known exactly |9. llCtlll| . Furthermore a theory of boundary induced phase transitions exists, which explains 
the phase diagram quantitatively in terms of the dynamics of shocks . In the stationary state these shocks exist 
as an upward density shock along the coexistence line between the high and low density phases, i.e., they connect a 
region with low density to the left of the shock position with a high density region to its right. The shock performs 
symmetric random walk between the boundaries of the system. 
. One may equip the system with the additional feature of local particle creation at empty sites with rate u a and 
annihilation with rate ujd (see Fig. ^) 8, 13]. In the thermodynamic limit L — > oo there are three regimes to be 
distinguished: If cu a and u>d are of an order larger than 1/L the steady state of the system will be that of Langmuir 
O ■ kinetics, i.e., there will be a uniform density of K = u> a /(u> a + LJd) m the system. In case of u> a and ujd being of smaller 
order than 1/L, the local kinetics is negligible and the system will behave as the TASEP. The case of the local rates 
being of order 1/L is the most interesting one and will be investigated further on. Writing 



uj a = tt a /L, LUd^^d/L (1) 

the phase diagram can be formulated in terms of f2 a , f2<z, p_ and p + . In Ref.^it was shown that for Q, a and Sl^ fixed, 
the phase diagram as a function of p_ and /?+ does not only exhibit the low-density and high-density phases known 
from the TASEP, but also a high-low coexistence phase. In this phase the shock does not move in the system but its 
position is a function of the rates p_ and p + (see Fig. 01. 

Parmeggiani et.al. presented not only Monte Carlo simulations but derived also a mean field equation for the density 
profile which was shown to coincide with the simulation profiles. We argue here that the mean field approximation 
can not be used in general. The coincidence with the Monte Carlo in Ref . |3 is due to lack of correlations in true steady 
state of the TASEP. We claim that the stationary density profile can be derived in general using a hydrodynamic 
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equation taking correlations into account (in case of the TASEP this equation is equal to that obtained with a mean 
field approach) . For the Katz-Lebowitz-Spohn (KLS) model, which is a generic model of interacting driven diffusive 
systems 0, we show that this hydrodynamic equation correctly describes the density profiles on a quantitative 
level, while a mean field approach would fail to reproduce even basic qualitative features of the system, e.g., phase 
separation into three distinct density regimes. 



II. HYDRODYNAMIC EQUATION 



In the following we are interested in the L — > oo limit which we perform by tuning the lattice spacing a = l/L — > 
and rescaling of time t = iiatticc/-^ (Eulerian scaling) to get the continuous (hydrodynamic) limit of the model. In 
this framework tt a ,d are the attachment/detachment rates per unit length. We claim that the hydrodynamic equation 
describing the time dependence of the local density p(x) for a general driven diffusive system with Langmuir kinetics 
takes the form 

-4,), e) 

where j(p) is the exact current in a driven diffusive system with homogeneous density p without LK and C(p) the 
source term describing the Langmuir kinetics. Here, we consider only that choice of C(p) which corresponds to the 
process depicted on Fig. ^ 

c(p) = n a (i-p(x,t))~n dP (x,t) (3) 

Other choices of C(fl), which might e.g. describe the local annihilation of particle pairs are to be discussed in a 
forthcoming publication [l6| . 

As is usually done in the rigorous derivation of the hydrodynamic limit of conservative systems |17| . our noncon- 
servative Eq. J3J implicitly assumes that the system is locally stationary because the exact form of the stationary flux 
is used. We argue that this assumption is justified since the nonconservative part of the dynamics of the system 
at macroscopic scale is so slow that locally the system reaches stationarity with respect to the conservative part of 
the dynamics. Any finite perturbation caused by the nonconservative dynamics would travel a macroscopic distance 
and hence dissipate before interacting with another perturbation. Hence the hydrodynamic description (after time 
rescaling t — > et) is adequate for describing the full dynamics. For physical insight in the formation of shocks one 
needs other tools which are discussed below. 

Rewriting equation J3J by using that dtp(x,t) — in the stationary state and d x j = dj/dp ■ dp/dx yields for the 
stationary density profile p(x): 

v c (p) d J^=C(p). (4) 

Here, v c = dj/dp is the collective velocity, i.e., the drift velocity of a center of mass of a local density perturbation 
on a homogeneous stationary background with the density p (for system with the Langmuir kinetics switched off) 
0,^2- The stationary density profile has to satisfy as well as the boundary conditions p(0) = p_ and p{X) = p+. 
As equation Q is of the first order there will be in general no smooth solution fitting both boundary conditions. In 
the original lattice model this discrepancy is resolved by appearance of shocks and/or boundary layers. To regularize 
the problem, one can add to and correspondingly to a vanishing viscosity term 

where v > is of order of l/L. This term makes the hydrodynamic equation second order and ensures a smooth 
solution fitting both boundary conditions. The shock has then a width of order l/L (see Ref. 0), i.e., in the thermo- 
dynamic limit the rescaled solution becomes discontinuous. We claim that equation 10) gives the same result in the 
L — > oo limit as the Monte Carlo, therefore it can be used as a tool to compute the stationary density profile. The 
main difference between (JSJ) and the MC is that the former does not take fluctuations into account which leads to a 
shock width of order l/L while in a MC after averaging it is of the order of 1/ \f~L due to the fluctuation of the shock 
position. 

The stationary density profile for a given j(p) and parameters Q a , f^, p_ and p + can be derived from the flow-field 
of the differential equation Q by using the rules, formulated and explained below: 
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(A) In the interior of the lattice the stationary density profile either follows a line of the flow field of the differential 
equation or makes a jump. Jumps can only occur between densities yielding the same current, i.e., the 
current is continuous in the interior of the lattice. 

(B) Let p'± be defined as limiting left and right densities with the boundary layers cut away: 

P-= lim p(x), p + = lim p(x), 

x^+Q x^l — 

where p(x) is the stationary profile in the hydrodynamic limit. The boundary layer at x = ( i.e., if p_ ^= p'_) 
has to satisfy the following condition: 

if p_ < p'_ then j(p) > j(p'_) for any p e (p_,p'_) (6) 
if p_ > p'_ then j(p) < j(p'_) for any p e (p'_,p_) (7) 

The condition for the stability of the boundary layer at x — 1 (if there is) is similar: 

if P+ < P+ then j(p' + ) < j(p) for any p e {p' + ,p+) (8) 
if p' + > p+ then j(p' + ) > j{p) for any p e (p+,p+) (9) 



(C) Shocks between a density p; to the left of the shock and p r to the right of the shock are stable only if they are 
stable in the absence of Langmuir kinetics pi Il8j. 



Remarks: 



Although LK does not conserve locally the number of particles, Eq. with the vanishing viscosity added (JSJ 
can be rewritten formally in the form 



at 



|-j(z, t) = 0, j(x, t) = j(p) - J C(p)dx - ^ - Ht) (io) 



where is some time-dependent function. Suppose that there is a shock at the position Xq connecting the 
densities pi and p r . The mass transfer across the shock is 

4 / ° + p(x,t)dx=j(X o +Q,t)-j(X o -0,t)=j(p r )-j(pi), (11) 

Ot Jx -0 

since the Langmuir term and the viscosity term change only infinitesimally across the shock. In the stationary 
state, the RHS of l|ll|l vanishes which explains the rule (A). 

• The rule (B) is due to the fact that in the boundary layer of the vanishing length SI — * 0, the LK term in ltTU|) 
can be neglected. Consequently, for the stationary current at the boundaries we have j(x) = j(p(x)) — = J, 
which yields the known maximization/minimization nrincinle[ll.l2lj|. and is equivalent to rule (B). Indeed at the 
left boundary J = j(p'J) (see (0 for notations), and if, e.g., p_ < p'_, then §f > 0. Consequently, we obtain 
j(p-) = J + v '^ > J, which is exactly @). Analogously one obtains Q-©. 

• The rule (C) is explained by the marginal role the Langmuir kinetics plays locally in space and in time. The first, 
LK is very slow locally for large L (see and the second, it acts "orthogonally" on the particle distribution, 
not affecting directly the particle motion. Hence, the local perturbations will still spread with the velocity 
corresponding to the local density level p, thus rendering the same stability conditions for a shock as for the 
diffusive system without LK. 

Condition (C) is easy to check geometrically through the current-density relation: an upward (downward) 
shock is stable if the straight line connecting the points (pi,j(pi)) and (p r ,j(pr)) stays below (above) the j(p) 
curve 0, . Because of criterion (A) these lines are always horizontal in this case which gives zero mean 
velocity (but not localization) for the shock in absence of Langmuir kinetics. 

• In the cases we have considered (ASEP, KLS model), the rules (A)-(C) define an unique stable solution (see an 
Appendix) and we believe that this is true also in general case, i.e., for arbitrary j(p) dependence and given 
choice © of Langmuir kinetics. 



In the following we apply the general theory to specific models. 
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III. REVISITING THE ASEP WITH LANGMUIR KINETICS 

Using the differential equation (0} and the rules given above we reconsider the ASEP with Langmuir kinetics [8Lll3|. 
Here, the current-density relation is given by j(p) — p(l — p), which yields v c {p) — 1 — 2p. Thus equation becomes 



(i - 2 P {x))d x p(x) = n„ - (n„ + n d ) P (x), (12) 

which is identical with the mean held equation in Ref. |8| in the thermodynamic limit. We would like to stress that 
this coincidence is caused by the fact that the mean held current-density relation for the TASEP is exact. As is 
demonstrated below, equation J3J also holds when this is not the case, as e.g. for the one-dimensional KLS model. 

Due to rule (A) as stated above (continuity of the current in the interior of the lattice) shocks in the interior can 
only occur in the case that pi = 1 — p r , as j(p) is symmetric to p = 1/2. Rule (C) (stability of the shock) furthermore 
requires that p r > pi- These observations coincide with the findings of 

We also applied our rules to k-hop exclusion models ( with LK added) , which are a generalization of the TASEP 
with stationary product measures and asymmetric current-density relations. Due to this fact shocks appear, which 
are non-symmetric with respect to p = 1/2. MC simulations are in full accord with our predictions |20). 



IV. KLS MODEL WITH LANGMUIR KINETICS 



A much studied one-dimensional driven diffusive system with interactions between the particles is the following 
variant of the KLS model H El EH : 

In the interior, particles at site i move to site i + 1, provided it is empty, with a rate that depends on the state of 
sites i — 1 and i + 2. 
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with 
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with 



rate 
rate 
rate 
rate 
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At site 1 particles can enter the lattice provided the target site is empty. The rate depends on the state of site 2. 
Similarly, particles can leave the system at site L with a rate depending on the state of site L—l. The boundaries mimic 
the action of reservoirs with densities p_ and p + . For p_ = p + the stationary state is that of an one-dimensional Ising 
model with boundary helds. The current-density relation can be calculated exactly using transfer matrix techniques 
|18| . It turns out that for strong enough repulsion between the particles (e > 0.9) a current-density relation with 
two maxima arises (see Fig- EJ • The parameter 8 determines the skewness of j(p) with respect to the vertical line 
p = 1/2. For 8 — 0, the system has particle-hole symmetry resulting in j(p) being symmetric with respect to 1/2. 
For simplicity we consider this case in the rest of the paper. 

The phase diagram of this family of models with strong particle repulsion is known to exhibit 7 different phases, 
among them two maximal-current and one minimal-current phase. The phase diagram is determined by the interplay 
of diffusion, branching and coalescence of shocks |21|. 

When equipping these models with Langmuir kinetics one expects that a very rich phase diagram with many more 
than the original 7 phases will appear. We will not attempt to give this full phase diagram here, but instead present 
two new features, which cannot be observed in systems without a concave region in the current-density relation: 
localized downward shocks and double shocks. 



A. Localized downward shocks 



In the regime where the current-density relation of the KLS model exhibits two maxima at densities p\ and p\, 
where p* < p% and a minimum at p = 1/2 (at 8 = 0) there is a region where downward shocks are stable according 
to Ref. [LSI I2ll 7and rule (C)). These are characterized by p\ <E (0.5, p£) and p r £ (p*, 0.5). This suggests that localized 
downwards shocks may appear when introducing the kinetic rates. In deed, in the KLS model with Langmuir kinetics 
for certain values of the boundary densities p_ and p+ , which strongly depend on the kinetic rates tt a and fi^, one 
gets a stable downward shock according to rules (A,B,C). We give an example for this case on Fig. 21 
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One can see that employing the general theory described above yields a stationary profile with a localized downward 
shock, which coincides with the MC results up to finite size effects, while a simple mean field approach would fail as 
it would not be able to capture the difference between the KLS model with e > and the TASEP (KLS with e = 0). 

B. Localized double shocks 

Let Pi 2 be defined as the inflection points of the current-density relation (p[ < p' 2 ). As is known from the studies 
of the KLS model if we start an infinite system from a step-like initial density profile with p_ e (pi, p' x ) on 

the left and p + S (p' 2 , p-i) on the right, we get a time-dependent solution having two shocks: One of these has negative 
mean velocity, while the other has positive and in the middle there is an expanding region with p = 1/2 (for 5 = 0) 
which corresponds to the minimal current phase in a system with open boundaries |l 81 12 1| . 

This leads us to the conjecture that introducing the kinetic rates for certain values of p_, p+, fi a , fid one may achieve 
a stable double shock structure. In Fig.[f)]we present an example for such a case. Application of rules (A,B,C), which 
is presented in detail in Appendix B yields the same double shock structure as the MC up to finite size effects. Note, 
that a simple mean field approach could not predict a double shock. 

V. CONCLUSIONS 

In this work we present a hydrodynamic equation which, together with some rules treating the discontinuities, 
correctly describes the stationary states of one-dimensional driven diffusive systems with Langmuir kinetics and 
open boundaries. It captures both systems without correlations in a steady state ( as e.g. the TASEP and the 
fc-hop exclusion models) and systems with correlations as the KLS model. For the latter the two new phenomena 
of a stationary localized downward shock and a localized double shock (corresponding to phase separation to three 
distinct regions) were presented which a mean field approach would not reproduce. The exact currect of driven 
diffusive systems without LK enters the hydrodynamic description since the bulk has sufficient time to relax between 
subsequent annihilation/creation events. An interesting, paradoxical feature of these phenomena is that fluctuating 
shocks get localized due to extra noise (LK), which is highly unexpected. 
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Appendix: DOUBLE SHOCK DENSITY PROFILE FROM THE RULES (A)- (C) 

Here we demonstrate how one determines the stationary density profile using the rules (A), (B) and (C) from the 
section |n] As an example we take the parameters which yield a double (localized) shock structure in the KLS model 
(p- = 0.23, p+ = 0.745, Q a = 0.03 and Q d = 0.01). The KLS-model parameters are: S = 0, e = 0.9 (see section ITvTl . 

First suppose that there is a boundary layer at x = 0. According to rule (B) it is stable only if p'_ > 1 — /?_ = 0.77. 
If this is the case then in the bulk there is no allowed jump since these trajectories of the flow-field (see Fig. [SJl stay 
always above p = 0.75 (rules (A) and (C)) which yields p' + > 0.75. But then the boundary layer at x = 1 does not 
satisfy rule (B). This contradiction shows that there is no boundary layer at x = 0. One can use the same argument 
to show that there is no boundary layer at x = 1 either. 

Now one can see that the stationary density profile close to the left boundary follows the line of the flow-field for 
which p(x = 0) = p- = 0.23. Since there is no boundary layer at the right end it is clear that somewhere in the bulk 
it has to make a jump. 

Note that this trajectory crosses the line p = p\ at x = x%. Suppose that the jump takes place before at x < x\. In 
this case, according to rule (A) it would jump over p% = 1 — p\ which would result in a boundary layer at x = 1 which 
is not allowed. If the jump takes place at x > x\ then p\ < p r < 0.5 and since from this region there is no allowed 
jump it would end up at pi < p', < 0.5 resulting again in an unstable boundary layer on the right side. This shows 
that the jump is located at x = x\ and from here the density profile follows the trajectory which starts at x = x\ 
with the value p = 0.5 + 0. 
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One can easily see that we need another jump to connect this trajectory with the one which ends at x = 1 with 
p = p + . Applying rule (A) (continuity of the current) we can get the point x 2 where the second jump is located. 
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Figure 1: Possible processes and their rates in the model of the ASEP with Langmuir kinetics 
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Figure 2: Plot of an average density of particles p versus rescaled coordinate x (site number /L) of a localized density shock 
in the ASEP with Langmuir kinetics. Parameters are p_ = 0.2, p+ — 0.6, Q a = 0.3 and Qd = 0.1. We show the results of both 
Monte Carlo simulations for L=100 




Figure 3: Current-density relation for the one-dimensional KLS model for various e. 
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Figure 4: Density of particles p versus rescaled coordinate x (site number / L) in a localized downward shock in the KLS model 
with Langmuir kinetics. Parameters are p_ = 0.64, p+ = 0.35, fl a = fid = 0.05. We show the results of both hydrodynamic 
equation and Monte Carlo simulation for L — 1000. The smoothness of the MC result is due to the fluctuation of the shock 
position^]. 




Figure 5: Path in the current-density relation for the profile shown in figure 2] 
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Figure 6: Density of particles p versus rescaled coordinate x (site number / L) in a localized double shock in the KLS model 
with Langmuir kinetics. Parameters are p_ = 0.23, p+ = 0.745, Q a = 0.03 and fid = 0.01. We show the results of both 
hydrodynamic equation and Monte Carlo simulation for L = 1000. The smoothness of the MC result is due to the fluctuation 
of the shock position|T^|. 




Figure 7: Path in the current-density relation for the profile shown in figure El 
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Figure 8: The flowfield of the hydrodynamic equation in the KLS model with Langmuir kinetics. Parameters are 5 — 0, e = 
0.9, Q a = 0.03, fid = 0.01. The thick lines show the stationary density profile for p_ = 0.23, p+ = 0.745 given by the rules 
(A,B,C). The dotted lines are p = pi ~ 0.24821, p = p2 ~ 0.75178 (see the subsection IIV Bl for notations). The axes: x is a 
rescaled coordinate [site number /L), p(x) is an average density of particles at point x. 



